{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<img src=\"../common/rfsoc_book_banner.jpg\" alt=\"University of Strathclyde\" align=\"left\">"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Notebook Set E\n",
    "\n",
    "---\n",
    "\n",
    "## 02 - QAM Modulation <a class=\"anchor\" id=\"QAM\"></a>\n",
    "\n",
    "Quadrature Amplitude Modulation (QAM) is a common form of modulation in wireless communication and can be seen in many modern radio communications systems. This scheme uses the principle of orthogonality to transmit two information signals within the same bandwidth. One signal is modulated on a sine wave at $f_{c}$ Hz, and a second signal is independently modulated on to a cosine wave, also at $f_{c}$ Hz (i.e. 90 degrees out of phase, or orthogonal). \n",
    "\n",
    "## Table of Contents\n",
    "* [1. Introduction](#introduction)\n",
    "* [2. QAM Modulation](#QAM_mod)\n",
    "* [3. QAM Demodulation](#QAM_demod)\n",
    "* [4. Phase Error With QAM](#QAM_phase_error)\n",
    "* [5. Conclusion](#conclusion)\n",
    "\n",
    "## Revision\n",
    "* **v1.0** | 22/12/22 | *First Revision*\n",
    "\n",
    "---\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 1. Introduction <a class=\"anchor\" id=\"introduction\"></a>\n",
    "Quadrature Amplitude Modulation can be used to achieve more bandwidth-efficient signalling. Using this approach, two signals are transmitted on a single complex carrier. The RF QAM transmitter can be seen in Figure 1.\n",
    "\n",
    "<figure>\n",
    "<img src='./images/QAM_modulation.svg' height='60%' width='60%'/>\n",
    "    <figcaption><b>Figure 1: QAM modulation architecture.</b></figcaption>\n",
    "</figure>\n",
    "\n",
    "For this notebook we will be using NumPy for computation and MatplotLib for Matlab-like visualisation of our waveforms. Let us begin by importing these libraries."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 2. QAM Modulation <a class=\"anchor\" id=\"QAM_mod\"></a>\n",
    "Firstly, we must establish the parameters that will be used throughout the notebook. These parameters match those seen in the amplitude modulation notebook, with the addition of a second baseband frequency for our second information signal."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Set basic params\n",
    "fs = 4096e6 # sample rate\n",
    "fb = 64e6 # frequency of baseband signal\n",
    "fb2 = 32e6 # frequency of second baseband signal\n",
    "A1 = 2 # baseband signal amplitude\n",
    "A2 = 1 # second baseband signal amplitude\n",
    "N_fft = 2048 # fft size\n",
    "\n",
    "t = np.arange(N_fft)/fs #time scale"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We now describe two information signals. For this demonstration, $g_{1}(t)$ is the same tone used previously in the amplitude modulation notebook, whilst $g_{2}(t)$ is a tone at a lower frequency. In this way, both tones can be transmitted using the same bandwidth as shown in the previous example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g1 = A1*np.cos(2*np.pi*fb*t)\n",
    "g2 = A2*np.cos(2*np.pi*fb2*t)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Calculate FFT\n",
    "g1_fft_result = np.fft.fft(g1, N_fft)\n",
    "g2_fft_result = np.fft.fft(g2, N_fft)\n",
    "\n",
    "# Get the corresponding frequencies, that depend on N_fft and Fs - freq. domain x axis\n",
    "freqs = np.fft.fftfreq(N_fft,1/fs)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These two signals can be plotted together on the time and frequency domains."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots(1,2, figsize=(15,4))\n",
    "axs[0].plot(g1[:200])\n",
    "axs[0].plot(g2[:200])\n",
    "axs[0].set_title('Time Domain Baseband Information Signals')\n",
    "axs[0].set_xlabel('Samples')\n",
    "axs[0].set_ylabel('Amplitude')\n",
    "axs[0].legend(('Information Signal 1', 'Information Signal 2'))\n",
    "\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(g1_fft_result[:int(N_fft/2)]))\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(g2_fft_result[:int(N_fft/2)]))\n",
    "axs[1].set_title('One Sided FFT plot of the Baseband Information Signals')\n",
    "axs[1].set_ylabel('Magnitude')\n",
    "axs[1].set_xlabel('Frequency, MHz')\n",
    "axs[1].legend(('Information Signal 1', 'Information Signal 2'))\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These two information signals are mixed (modulated) using In Phase and Quadrature carriers to give modulated signals that are 90 degrees out of phase.\n",
    "\n",
    "$$\n",
    "s_{1}(t) = g_{1}(t) \\cos(2 \\pi f_{c} t)\n",
    "$$\n",
    "\n",
    "$$\n",
    "s_{2}(t) = -g_{2}(t) \\sin(2 \\pi f_{c} t)\n",
    "$$\n",
    "\n",
    "The carrier frequency, $f_{c}$, is the same as that used in the previous example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fc = 1400e6 # Carrier Frequency \n",
    "\n",
    "# Find I and Q modulated signals\n",
    "s1 = g1*np.cos(2*np.pi*fc*t)\n",
    "s2 = -g2*np.sin(2*np.pi*fc*t)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find FFT of modulated signals\n",
    "s1_fft_result = np.fft.fft(s1, N_fft)\n",
    "s2_fft_result = np.fft.fft(s2, N_fft)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now plot the I and Q modulated signals. It is clear in the frequency domain that the bandwidth is $2f_{b}$, as we can see maximum and minimum frequency responses at $f_{c} - f_{b1}$ and $f_{c} + f_{b1}$, respectively."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots(1,2, figsize=(15,4))\n",
    "axs[0].plot(s1[:200])\n",
    "axs[0].plot(s2[:200])\n",
    "axs[0].set_title('Time Domain Modulated I and Q Signals')\n",
    "axs[0].set_xlabel('Samples')\n",
    "axs[0].set_ylabel('Amplitude')\n",
    "axs[0].legend(('I Modulated Signal', 'Q Modulated Signal'))\n",
    "\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(s1_fft_result[:int(N_fft/2)]))\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(s2_fft_result[:int(N_fft/2)]))\n",
    "axs[1].set_title('One Sided FFT plot of the Modulated I and Q Signals')\n",
    "axs[1].set_ylabel('Magnitude')\n",
    "axs[1].set_xlabel('Frequency, MHz')\n",
    "axs[1].legend(('I Modulated Signal', 'Q Modulated Signal'))\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These signals can then be summed together to create one signal, $y(t)$, which contains both information signals. \n",
    "\n",
    "$$\n",
    "y(t) = g_{1}(t) \\cos(2 \\pi f_{c} t) - g_{2}(t) \\sin(2 \\pi f_{c} t)\n",
    "$$\n",
    "\n",
    "In this way, two signals can be transmitted using a bandwidth of $2f_{b}$.\n",
    "\n",
    "Within Python, this is a simple case of summing together the two modulated signals, $s_{1}(t)$ and $s_{2}(t)$, as follows."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find sum of modulated signals\n",
    "y = s1 + s2 # sum of modulated signals"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find FFT of Summed Signal\n",
    "y_fft_result = np.fft.fft(y, N_fft)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This signal can now be plotted as before. We see that in both time and frequency domains the new signal $y(t)$ is the summation of $s_{1}(t)$ and $s_{2}(t)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots(1,2, figsize=(15,4))\n",
    "axs[0].plot(y[:200])\n",
    "axs[0].set_title('Time Domain Modulated Signal')\n",
    "axs[0].set_xlabel('Samples')\n",
    "axs[0].set_ylabel('Amplitude')\n",
    "\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(y_fft_result[:int(N_fft/2)]))\n",
    "axs[1].set_title('One Sided FFT plot of the Modulated Signal')\n",
    "axs[1].set_ylabel('Magnitude')\n",
    "axs[1].set_xlabel('Frequency, MHz')\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 3. QAM Demodulation <a class=\"anchor\" id=\"QAM_demod\"></a>\n",
    "\n",
    "At the receiver, quadrature demodulation involves multiplying the received signal by sine and cosine terms. The RF QAM receiver can be seen in Figure 2.\n",
    "\n",
    "<figure>\n",
    "<img src='./images/QAM_demodulation.svg' width='70%'/>\n",
    "    <figcaption><b>Figure 2: QAM demodulation architecture.</b></figcaption>\n",
    "</figure>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the In Phase channel, I, the output of the mixer is,\n",
    "\n",
    "$$\n",
    "x_{1}(t) = y(t) \\cos(2 \\pi f_{c} t),\n",
    "$$\n",
    "\n",
    "$$\n",
    "x_{1}(t) = 0.5 g_{1}(t) + 0.5g_{1}(t) \\cos(4 \\pi f_{c} t) - 0.5 g_{2}(t) \\sin (4 \\pi f_{c} t).\n",
    "$$\n",
    "\n",
    "Similarly, for the Quadrature Phase channel, Q, the output after the mixer is,\n",
    "\n",
    "$$\n",
    "x_{2}(t) = y(t)( -\\sin(2 \\pi f_{c} t)),\n",
    "$$\n",
    "\n",
    "$$\n",
    "x_{2}(t) = 0.5g_{2}(t) - 0.5g_{1}(t) \\sin(4 \\pi f_{c} t) + 0.5g_{2}(t) \\cos (4 \\pi f_{c} t).\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find demodulated signals - y(t)*(I or Q carrier)\n",
    "x1 = y*np.cos(2*np.pi*fc*t) # I demodulated signal\n",
    "x2 = y*(-np.sin(2*np.pi*fc*t)) # Q demodulated signal"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find FFT of demodulated signals\n",
    "x1_fft_result = np.fft.fft(x1, N_fft)\n",
    "x2_fft_result = np.fft.fft(x2, N_fft)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The resulting plots of the following code cell demonstrates that the desired baseband frequency components have been reconstructed. However, unwanted high frequency components are still present in the signal."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots(2,2, figsize=(15,4))\n",
    "axs[0, 0].plot(x1[:300])\n",
    "axs[0, 0].set_title('Time Domain Demodulated I Signal')\n",
    "axs[0, 0].set_xlabel('Samples')\n",
    "axs[0, 0].set_ylabel('Amplitude')\n",
    "\n",
    "axs[0, 1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(x1_fft_result[:int(N_fft/2)]))\n",
    "axs[0, 1].set_title('One Sided FFT plot of the Demodulated I Signal')\n",
    "axs[0, 1].set_ylabel('Magnitude')\n",
    "axs[0, 1].set_xlabel('Frequency, MHz')\n",
    "\n",
    "axs[1, 0].plot(x2[:300])\n",
    "axs[1, 0].set_title('Time Domain Demodulated Q Signal')\n",
    "axs[1, 0].set_xlabel('Samples')\n",
    "axs[1, 0].set_ylabel('Amplitude')\n",
    "\n",
    "axs[1, 1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(x2_fft_result[:int(N_fft/2)]))\n",
    "axs[1, 1].set_title('One Sided FFT plot of the Demodulated Q Signal')\n",
    "axs[1, 1].set_ylabel('Magnitude')\n",
    "axs[1, 1].set_xlabel('Frequency, MHz')\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Low pass filters are used to remove the high frequency terms so that the information can be perfectly recovered, as with amplitude modulation. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Define lowpass filter\n",
    "f_cutoff = 0.1 # Cutoff frequency as a fraction of the sampling rate\n",
    "b = 0.08  # Transition band, as a fraction of the sampling rate (in (0, 0.5)).\n",
    "\n",
    "N = int(np.ceil((4 / b)))\n",
    "if not N % 2: N += 1  # N is odd.\n",
    "n = np.arange(N)\n",
    "\n",
    "h = np.sinc(2 * f_cutoff * (n - (N - 1) / 2)) # Compute sinc filter.\n",
    "w = np.blackman(N) # Compute Blackman window.\n",
    "h = h * w # Multiply sinc filter by window.\n",
    "h = h / np.sum(h) # Normalize to get unity gain."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Apply filter to demodulated signals using convolution\n",
    "z1 = np.convolve(x1, h)\n",
    "z2 = np.convolve(x2, h)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find FFT of filtered signals\n",
    "z1_fft_result = np.fft.fft(z1, N_fft)\n",
    "z2_fft_result = np.fft.fft(z2, N_fft)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can now plot the filtered signals in the time and frequency domains. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots(1,2, figsize=(15,4))\n",
    "axs[0].plot(z1[:200])\n",
    "axs[0].plot(z2[:200])\n",
    "axs[0].set_title('Time Domain Lowpass Filtered Demodulated I and Q Signals')\n",
    "axs[0].set_xlabel('Samples')\n",
    "axs[0].set_ylabel('Amplitude')\n",
    "axs[0].legend(('I Demodulated Signal', 'Q Demodulated Signal'))\n",
    "\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(z1_fft_result[:int(N_fft/2)]))\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(z2_fft_result[:int(N_fft/2)]))\n",
    "axs[1].set_title('One Sided FFT plot of the Lowpass Filtered Demodulated I and Q Signals')\n",
    "axs[1].set_ylabel('Magnitude')\n",
    "axs[1].set_xlabel('Frequency, MHz')\n",
    "axs[1].legend(('I Demodulated Signal', 'Q Demodulated Signal'))\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Notably, the input information signals have both been preserved. However, the amplitude of the signals have halved, which follows is the same property we observed in the previous example. The amplitude is halved because QAM is still a form of amplitude modulation and we accept that there will be an amplitude loss upon reconstruction. The amplitude loss can be easily corrected by applying a gain.\n",
    "\n",
    "We can see that the final signals can be represented as,\n",
    "\n",
    "$$\n",
    "z_{1}(t) = 0.5g_{1}(t),\n",
    "$$\n",
    "\n",
    "$$\n",
    "z_{2}(t) = 0.5g_{2}(t).\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4. Phase Error With QAM <a class=\"anchor\" id=\"QAM_phase_error\"></a>\n",
    "\n",
    "If the receiver's local oscillator is shifted by $\\theta$, with respect to the received signal, then the I and Q output signals will be mixed. \n",
    "\n",
    "For the I channel, the phase shift is defined as,\n",
    "\n",
    "$$\n",
    "x_{1}(t) = y(t) \\cos(2 \\pi f_{c} t + \\theta),\n",
    "$$\n",
    "\n",
    "$$\n",
    "x_{1}(t) = 0.5 [g_{1}t \\cos(\\theta) + g_{2}(t) \\sin(\\theta)].\n",
    "$$\n",
    "\n",
    "For the Q channel, the phase shift is defined as,\n",
    "\n",
    "$$\n",
    "x_{2}(t) = y(t) (-\\sin(2 \\pi f_{c} t + \\theta)),\n",
    "$$\n",
    "\n",
    "$$\n",
    "x_{2}(t) = 0.5 [-g_{1}t \\sin(\\theta) + g_{2}(t) \\cos(\\theta)].\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create phase shift\n",
    "phaseshift = (np.pi/3) #60 degree phase shift\n",
    "\n",
    "# Apply phaseshifted I and Q carriers to modulated signal y(t)\n",
    "x1_phaseshift = y*np.cos(2*np.pi*fc*t + phaseshift)\n",
    "x2_phaseshift = y*(-np.sin(2*np.pi*fc*t + phaseshift))\n",
    "\n",
    "# Apply lowpass filter to I and Q demodulated phase shifted signals\n",
    "z1_phaseshift = np.convolve(x1_phaseshift, h)\n",
    "z2_phaseshift = np.convolve(x2_phaseshift, h)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Find FFT of phase shifted signal\n",
    "z1_phaseshift_fft_result = np.fft.fft(z1_phaseshift, N_fft)\n",
    "z2_phaseshift_fft_result = np.fft.fft(z2_phaseshift, N_fft)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From the plots below, it is clear that frequency components from $g_{1}(t)$ and $g_{2}(t)$ are present in both I and Q channels. In the time domain, it can be seen that the demodulated I and Q phases are no longer simple sine waves – in fact, each is now composed of two different frequencies. This is due to the phase shift in the receiver local oscillator causing the I and Q phases to be mixed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "fig, axs = plt.subplots(1,3, figsize=(20,4))\n",
    "axs[0].plot(z1_phaseshift[:200])\n",
    "axs[0].plot(z2_phaseshift[:200])\n",
    "axs[0].set_title('Time Domain Lowpass Filtered \\nDemodulated Signals with Phase Offset')\n",
    "axs[0].set_xlabel('Samples')\n",
    "axs[0].set_ylabel('Amplitude')\n",
    "axs[0].legend(('I Demodulated Signal', 'Q Demodulated Signal'))\n",
    "\n",
    "axs[1].plot(freqs[:int(N_fft/2)]/1e6, np.abs(z1_phaseshift_fft_result[:int(N_fft/2)]))\n",
    "axs[1].set_title('One Sided FFT plot of the I Lowpass Filtered \\nDemodulated Signals with Phase Offset')\n",
    "axs[1].set_ylabel('Magnitude')\n",
    "axs[1].set_xlabel('Frequency, MHz')\n",
    "\n",
    "axs[2].plot(freqs[:int(N_fft/2)]/1e6, np.abs(z2_phaseshift_fft_result[:int(N_fft/2)]))\n",
    "axs[2].set_title('One Sided FFT plot of the Q Lowpass Filtered \\nDemodulated Signals with Phase Offset')\n",
    "axs[2].set_ylabel('Magnitude')\n",
    "axs[2].set_xlabel('Frequency, MHz')\n",
    "\n",
    "fig.tight_layout()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The I and Q signals are interfering with each other due to the phase error. DSP receivers can be designed to calculate these phase errors and correct for them."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5. Conclusion <a class=\"anchor\" id=\"conclusion\"></a>\n",
    "\n",
    "Throughout this notebook, we have explored the QAM modulation scheme and investigated how it can be used to transmit two information signals using two orthogonal carriers at the same frequency. This technique improves the bandwidth efficiency of a radio communications system. In the next notebook, we will explore a complex representation of the QAM scheme and demonstrate its equivalence to the trigonometric version introduced here."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "---\n",
    "\n",
    "[⬅️ Previous Notebook](01_amplitude_modulation.ipynb) || [Next Notebook 🚀](03_complex_qam.ipynb)\n",
    "\n",
    "Copyright © 2023 Strathclyde Academic Media\n",
    "\n",
    "---\n",
    "---"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
